{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Elevation profiles extracted from SRTM over Mt Baker compared with ICESat-2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import geopandas as gpd\n",
    "import rasterio\n",
    "import matplotlib.pyplot as plt\n",
    "import glob\n",
    "from topolib import gda_lib\n",
    "import topolib"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Read reference DEM"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "dem_fn = '/home/jovyan/data/srtm_elevation/SRTM3/cache/srtm_wa_subset.vrt'\n",
    "dem = rasterio.open(dem_fn)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Read ATL06 data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_dir = '/home/jovyan/data/nsidc/**/'\n",
    "ATL06_list = sorted(glob.glob(data_dir + \"*.h5\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sample DEM at ICESat-2 ATL06"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "ATL06_fn = ATL06_list[0]\n",
    "dataset_dict={'land_ice_segments':['h_li',\n",
    "                                   'delta_time',\n",
    "                                   'longitude',\n",
    "                                   'latitude'],\n",
    "              'land_ice_segments/ground_track':['x_atc']}\n",
    "\n",
    "ATL06_gdf = gda_lib.ATL06_2_gdf(ATL06_fn,dataset_dict)\n",
    "ATL06_gdf = ATL06_gdf.to_crs(dem.crs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "poly = gda_lib.dem2polygon(dem_fn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "glas_gdf_aea_rgi = gpd.sjoin(ATL06_gdf, poly, op='intersects', how='inner')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "points_xy = list(zip(glas_gdf_aea_rgi.geometry.x, glas_gdf_aea_rgi.geometry.y))\n",
    "\n",
    "refdem_sample = []\n",
    "for val in dem.sample(points_xy):\n",
    "    refdem_sample.append(val[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "glas_gdf_aea_rgi['srtm_height'] = refdem_sample"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "glas_gdf_aea_rgi['diff'] = glas_gdf_aea_rgi['srtm_height'] - glas_gdf_aea_rgi['h_li']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD4CAYAAAAAczaOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO2df7CV1XnvP99zOJAj0SIRvIgQCBdNxR9Yz9XTC1JrS7GmUZOOUQpRr5lYHdsboyEJgYmaCY0NijaZMRmpTpqJATXq1jRpKLc3tA0DWnAjR1SqRitsuWqiThglyo/n/vG+G182+5y93v3r7Hfv5zNz5rxnvWvvvbau98taz3p+yMxwHMdxOoOu4R6A4ziO0zxc9B3HcToIF33HcZwOwkXfcRyng3DRdxzH6SBGDPcAKnHMMcfYlClThnsYTpuyefPmX5nZuGZ/rs9rp5EMNa9bXvSnTJnCpk2bhnsYTpsi6b+G43N9XjuNZKh57eYdx3GcDsJF33Ecp4Nw0Xc6lh07dgCcIOkZSdskfQ5A0kxJGyVtkbRJ0pnF10haLOl5SdslzUu0nyFpIL73LUlq/jdynMq46Dsdy4gRIwB2mtnvAv3AtZJOAr4J3GxmM4Gvxn8T37sUmAGcB9wpqTt+u+8AVwHT45/zmvhVHCcYF32nY5kwYQLAOwBmtht4BpgIGHBU3O13gFfi6wuB1Wb2rpm9CDwPnClpAnCUmW2wKJnV94GLmvZFHCcFLe+94zhDsTQ3wKrHdrDfjG6J+WdN4usXnZL6fSRNAU4HHgOuA9ZIupVoYfQ/424TgY2Jl+2M2/bG16XtpZ9xFdFugMmTJ6ceo9NZ5PIFbnp0G2/t2QvA0Uf0cOPHZ3DR6YdNrVT4St/JLEtzA/xg48vsjzPF7jfjBxtfZsHKDaneR9IHgQeB68zsN8A1wOfNbBLweeDuYtcyL7ch2g9tMLvLzPrMrG/cuKaHBjgZIpcvsOiBJw8KPsCb7+xl0Y+eJJcv1PTeFUVf0iRJPy9z2HVffNC1RdJLkrbE7VMk7Unc+27ivfywy6kbP3zs5bLt6194g6W5gdC3EZHg32tmD8VtlwPF6weA4kHuTmBS4rXHE5l+dsbXpe2OUxXL12xn74HD097v3W8sX7O9pvcOWenvA24oPewys0vMbGZ82PUg7z8kAC8U75nZ1Yl2P+xy6kaZZ+IgP9j4csUVUVxL4sPAM2a2InHrFeAP4utzgefi60eBSyWNkjSVaA4/bma7gN2S+uOFzGXAI6m/kOPEvPLWnqruhVBR9M1sl5k9EV8nD7sAiCf5p4BVQ72PH3Y59SRki1tpRbR+/XqADwHnJnam5wOfBW6T9CTwN8R2eDPbBtwPPA38DLjWzPbHb3cN8PdEh7svAP+U/ls5nU4uX2D6V35yuG0wwXFjemv6jFQHuSWHXUXOBl41s+cSbVMl5YHfAEvN7N+J/qGoeNgVf44feDlDcvOPt1XsU2lFNHv2bIDNZtZX5vYZ5V5jZsuAZWXaNwEnVxyU4wxCLl/guvu2DNmnS7Bo3ok1fU7wQW6Zw64i8zl0lb8LmGxmpwPXAz+UdBSBh13gB15OZd58Z2/FPiPcTcHJECG2+gMG1923pabD3KDHQlIPhx92IWkE8EngvmJb7MP86/h6M9FW9wT8sMupE6ETfu+B8L6OM9yksdXXIvwh3jsiclkrPewC+GPgWTPbmeg/rhilKOkjRIddv/TDLqdehJh2iix+aGsDR+I49SOtrb5aL56Qlf4s4NMcftgFUUh66QHuHGBrfAj2I+BqM3sjvueHXU7NhJh2iuzZe8BX+04mSGurr9aLp+JBrpn9gvL2eMzsijJtDxKZgsr198Mup+ksemBLzVGMjtNoLjp9Ipv+6w1+sLF8/EkpRmS+TDu3/ajLyRxjentS9XfbvpMFcvlCsOAXueGB9BG6LvpO5rjpghnpX/No+DmA4wwH1djo9x9IH6Hrou90BMkcJo7TilRro0/7Ohd9J3Ok8d5xnKxQS6RtGhOPi76TOYby3ukaJIXfyG7P7ee0NkeMrE6OjXS2fRd9p60YTNvf229pMm86TlNZsHIDz732dtWvT2Pbd9F3MsVQqxkReeoMRkjmTccZDta/8EblThUIte276DuZYsnDg6/Wh8pMWMS9eJx2JfRMwEXfyQy5fIG339s/6P2JY3oHtekXcS8eJ2tIYUL9hx8NS07pou9khi89OHQenUXzTuQvzvJU3E57YQYrLplZsd/Pn3096P1c9J1MkMsXeHff4Ab73p4uLjp9YlVF0R2nlZk4pjco1YLb9J22opJv/jc+eerB64mBts0dO3YAnFBa/xlA0l9L2h63fzPRvjiu8bxd0rxEu9d/dqpm1rSxZduTRVO6K0wpt+k7bUWlzJrJldBQ2QqTD86IESMAdpbWf5b0h8CFwKlmNgO4FUDSSUSZZWcQ1Xe+s5hGHK//7FTJ0txAWe+dUSO6WPGpmQfn9vyzJg36Ht1dCs7S6aLvZJ6F/Yfa8S86feKgK6fkgzNhwgSAd+Cw+s/XALeY2bvxvdfil1wIrI4LBb1IlCL8TK//7FTL0tzAoEnWJo891Kzz9YtOOWyuA4we2c1tF58WnG0zVY1cxxkuxvT2lPW86emirB3/3s/+PktzA6x6bAf7zeiWmH/WpEFt/iX1n5cDZ0taBvwW+IKZ/QfRPwgbEy8r1nneS0D9Z6/97JSy6rEdg94rF6z19YtOqfncykXfyQQ3XTCDRQ88yd4D73vj93SJ5RefNuhrQh+Q0vrPcRnQo4lMPv8DuD+uAjdYneeg+s9mdhdwF0BfX19IWIHT5uy35k8DF30nExS3rsvXbOeVt/Zw3JheFs07sR7FUcTh9Z93Ag/FpprHJR0Ajonbk4bVYp1nr//sVEW31HThd9F3MsNFp0+sawWsSNP5MPDvJfWfc8C5wDpJJwAjgV8BjwI/lLQCOI7owPZxM9svabekfiLz0GXAt+s2UKdt+ci4IwbNuTN9/OiGfKaLvpMJ0tjnQ1m/fj3Ah4jrP8fNXwHuAe6R9BTwHnB5vOrfJul+4GlgH3CtmRVDhK8Bvgf0EtV+9vrPzpAszQ0MKvhHjepm7fXnNORzXfSdlmfByg2HuLTtNzvo8VCL8M+ePRtgs5n1lbm9sNxrzGwZsKxMu9d/dlIx1CHu2+8NkTmwRtxl02lpcvnCoBkIh3poHKeVyeULQ9ryG2nnryj6kiZJ+nlp1KKkmyQVJG2Jf85PvMajFp26MFQk7nB4PjhOreTyBRY/NHRth0YKY4h5Zx9wg5k9IelIYLOktfG9283s1mTnkqjF44D/I+mE2PZZjFrcCPyUKGrRbZ/OoFSKxHWcrLF8zXb27B08WyyEpQmvloqib2a7gF3x9W5JxajFwTgYtQi8KKkYtfgScdQigKRi1KKLvnOQXL7A8jXbKVRZJNpxWpmluYFhn9upDnJLohZnAX8l6TJgE9Fu4E1qjFqMP8cjFzuQ0gPbSoQmVnOc4SSXL7DogS1DVnUrpVJytVoIPsgtjVokMtVMA2YS7QRuK3Yt8/LgqEWIIhfNrM/M+saNCysM4GSXpbkBpnz5J6lLxoUmmHKc4SKXL3DdfekEH4ZOrlYrQaIvqYeSqEUze9XM9pvZAWAlcGbc3aMWnWCGSjg1FAv7J9c1UMtx6k1R8NPS29PV0LoQFc07sYfN3cAzyahFSRNiez/AJ4Cn4muPWnQqUq3YF/FiKU6rUo05J8lvq31hICE2/VnAp4GBkqjF+ZJmEploXgL+EsDMPGrRGZJaBd9xWpVcvsD1922hFtkOLYZSLSHeO7+gvD3+p0O8xqMWnUG5t0bB9wNcp9VYmhvg3o0v1+xqmaYYSrV4Ggan6dTyYPR0N/6hcJw01HPnmqYYSrW46DtNZe6KdVW/9ugjerjx4zP8ANdpKeqVDiS0AHqtuOg7TWOorIKVuOOSmS72TstRKYdOKM0w6xTxhGtOU6hlCzxr2lgXfKflqNYlsxTRHLNOEV/pOw1n7op1Va3wBSzon+zumU5LsnzN9prfo7enm2988pSmLmp8pe80lGoEX0TmnBdv+VhDBX/Hjh0AJ5RmkD04DukLkkzSMYk2zyDrAFSVQ2fWtLFMHNOLiGz4zRZ88JW+00CqseGPHtnNsk8050EYMWIEwE4zOymZQdbMnpY0CZgLHLRJeQZZp8iClRtSv2Zhi+xaXfSdhpHGq2FMbw9bbvyTBo7mcCZMmADwDhyWQfZp4Hbgi8AjiZd4BlkHIFWeqFYR+yJu3nEaRhqvhpsumNHAkVQmmUFW0gVAwcyeLOk2EUj+S1bMFDuRgAyykq6StEnSptdff72Oo3dalVnTxraU4IOLvtNAulJYtYfTOyeZQZYodcgS4KvlupZpC84g69lj24NcvhDUr6cL7v3s7zd4NOlx0XcaRneg6C/sH9aaCeLQDLLTgKnAk7HZ5njgCUn/Dc8g2/EUc+uEsK+xedOqxkXfaQi5fCEoy+Bw2jstMj99mEQGWTMbMLPxZjbFzKYQCfrvmdn/I8oge6mkUZKm8n4G2V3Abkn9sdfOZRx6FuC0CcvXbA9OptboxGnV4ge5TkNY8vDQhZ9nTRs77Fvf9evXA3wIODeZQdbMyiYT9AyyTqibZpdat8iPi75TN0Lr2/b2dA274APMnj0bYLOZ9Q3WJ17tJ//2DLIdSJpSniO6xK1NjLBNi4u+UzO5fIHFD21lT2Dxh0YXiXCcepJG8LOQI8pF36mJXL7A5+/bkipdcqvaOh2nSDX58bullhd8cNF3auQLDzyZ6sHoaWI2QcephmqTAzaymHk9ce8dp2oWrNzAvgPp0soub2Fbp+NAdfnxP9CtlgvCGgwXfacqcvlCqlB0SBes5TjDRTX58X+7v/ac+s3CRd+pipt/vC31aw5Yda9znGaxNDe0q3E74KLvVMWb7+xt6uscp9HUs9ZtK1NR9CVNkvTz0pzjkpZLelbSVkkPSxoTt0+RtEfSlvjnu4n38pzjbUAtdW4dpxWpVfBnTRtbx9E0lpCV/j7gBjP7XaAfuDbOK74WONnMTgX+E1iceM0LZjYz/rk60V7MOT49/jmvHl/CaR4LVm6ous4tRCmUHaeVqIfgt0KwYSgVRd/MdpnZE/H1buAZYKKZ/bOZ7Yu7beTQhFOHIWkCcc5xi5KeFHOOOxlhaW4g9eFtEjH8KZQdp5RqvHWKLOyfnCnBh5Q2/WTO8ZJbV3JorpGpkvKS/lXS2XFbUM7x+HM873iLUQ975+0ZiFZ0OotcvlCVtw60XnGUUIKDs5I5x83sN4n2JUQmoHvjpl3AZDP7taQzgJykGQTmHIco7zhwF0BfX192fKHalFy+UJPg93TJ/fOdliOXL3BdYJrkJFkz55QStNKX1MOhOceL7ZcDfwYsiE02mNm7Zvbr+Hoz8AJwAp5zPJOkfTCmjx/NHZfMPKT4swu+04qkdR8WUW6dLAs+BKz0Yw+bu0nkHI/bzwO+BPyBmb2TaB8HvGFm+yV9hOjA9pdm9oak3ZL6icxDlwHfru/XcepJWsE/9siRrL3+HGB4K2E5TiVy+UIq9+HuLnFbmyxeQsw7s4BPAwPJnOPAt4BRwNrY83Jj7KkzB/iapH3AfuBqMyue/nnO8QxRKSd+kqzaN53OI+1iZtSILv72z09tC8GHANE3s19Q3h4/WKGJB4lMQeXuec7xDPH2e/srdyJyw8yi4O/YsQPgBEnPAAeAu8zs7yQtBz4OvEdknvxfZvYWgKTFwGeIFjT/28zWxO1n8P6C5qfA54omT6e1+NKDW4P7ZiFVclo8ItcpS2jxZ8iuG+aIESMAdobGoMT3LgVmEMWY3CmpO347j0HJAAtWbuDdwOK1C/snt53gg4u+U4ZcvsCiHz0Z1HfWtLGZfTAmTJgA8A4Ex6BcCKyOnRVeBJ4HzvQYlOwQGmfSzuZKF33nMJY8PMDegKyBWXddSxIYgzIRSEbyFGNNgmJQPP4kG2TVXBmKF1FxDiGXLwTZ8tvJ1pkiBmWwWJOgGBSPPxleFqzcENQvq+bKUFz0nUP4wgNhZp12EXwiwR4qBuWPEgeyO4FkeaRirInHoLQw7VbjtlbcvOMcJLQSVrskTYu1/MMMHoNyQTIGBXgUuFTSKElTiQ5sHzezXcBuSf1xXMtlwCPN+h7O4KQR/KzUuK0VX+k7Bwl9ONpl+7t+/XqADwHnhsSgmNk2SfcDTxOZfa41s6ItzGNQWpA0CQKzUuO2Vlz0nVS0kxvb7NmzATabWV/JrbIxKABmtgxYVqbdY1AyTDt765Tiou8AYX75nWDvdDqTThF8cJu+QyT4ix+qnHLBBd/JGiEVrbo7rICfi77D8jXb2bM3LOWC42SF0KI/nWLLL+Ki7/DKW3sq9mkXjx2nMwgt+tNJtvwiLvoOR4zsrtinXTx2nM4gtARipwk+uOh3PCERuL09XW7PdzJFSAnEzrLkv4+LfocTUj3oG588tQkjcZzmsqB/8nAPYVhw0e9wKlUPynIWTccZjN6ero407YCLfkezNFfZTbNdsmg6nUPIvP7t3rCc+u2Ii34Hc28F7wb32HGyRqjXznFjepswmtbERb9DyeULh+f+LcE9dpysUWkhA1GR80XzTmzCaFoTF/0O5SsPDV0n1G35ThaptJAZPbKb2y4+raPntufe6VDeqWDTdFu+kzVCiqRs+5qXLq640pc0SdLPJT0jaZukz8XtYyWtlfRc/PvoxGsWS3pe0nZJ8xLtZ0gaiO99K8497jSZuSvWDXn/6CPclu9kizR58zudEPPOPuAGM/tdoB+4VtJJwJeBfzGz6cC/xH8T37sUmAGcB9wpqRjy+R3gKqLiE9Pj+04TWbByA8+99vaQfW78uNvynewQmmMnJPlaJ1BR9M1sl5k9EV/vBp4hKvp8IfAPcbd/AC6Kry8EVpvZu2b2IvA8cKakCcBRZrYhLj/3/cRrnCYQuhrqFHvnjh07AE7wXWx2CfXWmTVtrJssY1Id5EqaApwOPAYcG5eJI/49Pu42EUgmvtgZt02Mr0vby33OVZI2Sdr0+uuvpxmiMwihgt9Jq6ERI0YA7PRdbDYJFXzwM6okwQe5kj5IVED6OjP7zRALmXI3bIj2wxvN7gLuAujr66ucRMMZklNv/Bm/eTcsdXInPRwTJkwAeAeiXayk5C72nLjbPwDriGrmHtzFAi9KKu5iXyLexQJIKu5ivWRiA8jlC1x335bKHWOOGlU5oWAnEbTSl9RDJPj3mtlDcfOrscmG+PdrcftOIJmg+njglbj9+DLtToPI5QtM+fJPggV/YYfmIoHm7WKd2kgr+ABbb/ZNV5IQ7x0BdwPPmNmKxK1Hgcvj68uBRxLtl0oaJWkq0Vb38fjh2S2pP37PyxKvcerMgpUbUj8cnZqLpHQXO1TXMm3Bu1g3W9ZGLl/g8ynn9B2XzGzQaLJLyEp/FvBp4FxJW+Kf84FbgLmSngPmxn9jZtuA+4GngZ8B15pZcal5DfD3RIe7L+Db34ZQjftaBz8cokm7WDO7y8z6zKxv3Lhx9f0WbU5xhZ/G1us1nctT0aZvZr9g8NTTfzTIa5YBy8q0bwJOTjNAJx25fCG14C/sn9yRD0fkRMaHgX8fZBd7C4fvYn8oaQVwHO/vYvdL2i2pn8g8dBnw7eZ8i84gJAV4kk6d0yF4RG6bsXzN9lT9O3k1tH79eoAPEe9i4+avEIn9/ZI+A7wMXAzRLlZScRe7j8N3sd8Deol2sL6LrRO5fKFiCvAknVgCMQ0u+m1GIaDeLcD08aNZe/05jR1MizN79myAzWbWV+a272JbgLQHt528iAnFRb+NyOULQf38wXCyQhqzjs/rMFz024gv/ujJin1euuVjTRiJ49SHULOO2/DD8dTKbcR7+4f2behgDx2njXEbfjpc9NuESiXifCXkZI1Qc6ULfjpc9NuESjlI/MFwskQuX+D6gAPcTo4irxYX/Tag0iq/k5KoOe3B8jXbqVS63M061eGi3wZUWuV3UhI1pz14JcD12AW/Olz0M06lVb7jZJHuLi9H0Chc9DNMLl+ouMp3m6eTNeauWMe+A55RvVG46GeYJQ9XXuX7FtjJGpXKeYIvZmrBRT/DvP3e0Hnye/z/rpMx5q5YV7GPH+DWhstCRvnvi39Ssc/yiz0Yy8kOc1esC1rlu+DXhqdhyCBnLVvLvgomz56uzilw7mSfs5at5dXd71Xsd+yRI5swmvbGV/oZY+6KdUEPh6/ynayQRvAfWzK3CSNqb3ylnyFCH47p40f7Kt/JDCFzGnDBrxO+0s8IC1ZuCHo4PtCtjs+T72SH0Pw6R43qbvBIOgcX/QyQpgTis8vOb/Bo2ocrr7wS4DRJTxXbJM2UtDGuBb1J0pmJe4slPS9pu6R5ifYzJA3E974lySOLAggtdH7UqG623nxeE0bUGbjoZ4DQykHuu5yOK664AuC5kuZvAjeb2Uzgq/HfSDoJuBSYAZwH3CmpuPz8DnAVUc3c6fF9ZwiW5gaCCp1PHz/aBb/OuOi3OCGumRAdcrkrWzrmzJkDUa3bJAYcFV//DvBKfH0hsNrM3jWzF4HngTMlTQCOMrMNFlVa/z5wUcMHn2GW5gYqRpJDJE5uqqw/FUVf0j2SXivZAt8Xb3+3SHqpWFRa0hRJexL3vpt4jW+BUxLimglRFk0/5Kob1wHLJe0AbgUWx+0TgR2JfjvjtonxdWn7YUi6KjYZbXr99dfrPvCssOqxHZU7QcVdgFMdISv971GyXTWzS8xsZrwFfhB4KHH7heI9M7s60e5b4BSEumb2dHkWzTpzDfB5M5sEfB64O24vt0ixIdoPbzS7y8z6zKxv3LhxdRls1sjlC+y3MDk/bkxvg0fTmVQUfTP7N6DsKWK8Wv8UsGqo9/AtcDpy+UJQZCK4P34DuJz3FzEPAMWD3J3ApES/44lMPzvj69J2p4TQwihFFs07sYGj6VxqtemfDbxqZsnDsKmS8pL+VdLZcVvwFhh8G7x8zfagfndcMtP98evPK8AfxNfn8v5B76PApZJGSZpKtFt93Mx2Absl9ceLoMuAR5o96CwQUhiliM/txlFrcNZ8Dl3l7wImm9mvJZ0B5CTNIMUWGKJtMHAXQF9fX8eZ9goBBSReuuVjTRhJezN//nyAjxJtWncCNwKfBf5O0gjgt0QmScxsm6T7gaeJDn+vNbNixrtriMygvcA/xT9Ogly+EDSvASaO6XXBbyBVi378UHwSOKPYZmbvAu/G15slvQCcgG+Bg1mwckPFPhPd1lkXVq1axerVq7eaWV/JrTPK9TezZcCyMu2bgJMbMMS2II1Zp7tLbtZpMLWs9P8YeNbMDpptJI0D3jCz/ZI+QrQF/qWZvSFpt6R+4DGiLfC3axl4O5DLF1j0wBb2hu55Y/yhcLJCMQArdLt+28Wn+Sq/wVQUfUmrgHOAY4pbYDO7myhQpfQAdw7wNUn7gP3A1WZWPAT2LXCCBSs3BEfZluIPhZMFiiv8UMHvlnxuN4GKom9m8wdpv6JM24NELpzl+vsWOGZpbqBqwR/T21Pn0ThOY/jSg1uDD24B5p81qXInp2Y8y+YwcG9ANGI5BNx0wYz6DsZx6kw1u1ivhtU8XPSbTC5fqDrS8HZ3Y3NamGpNlt2SC34T8dw7TSakmHk5FvZPdsF3WpZazqjcrNNcXPSbSC5fqFjMvByzpo31lZDT0lQr+G7WaT5u3mkiN/94W6r+wk06Tvvigj88+Eq/ibz5zt5U/V3wnSywNJfeZOmCP3y46DeJs5atTdV/1rSxLvhOyxOaGz+JH9wOLy76TSC0oHkST5fsZIHQ3PhJ/OB2eHHRbzC5fCG14HvZQycLzF2xLjg3fhE36ww/LvoNJrS+bZHp40f7Q+G0PHNXrAuu+VDEzTqtgYt+A/nokp+m6j9r2livCepkgrSCD27WaRXcZbNBLFi5gd/uD9v6jhrRxd/++al+cOtkgpD036W4Wad18JV+gwgJVhFRhaDtX/9TF/xh4MorrwQ4TdJTyXZJfy1pu6Rtkr6ZaF8s6fn43rxE+xmSBuJ734oraLUlnlcn+7joN4AQs86xR47kxVs+5mI/jFxxxRXwfjlEACT9IXAhcKqZzQBujdtPIkonPgM4D7hTUnf8su8QVdiaHv+c14ThDwtpBd+jyVsPN+/Umbkr1lU063ygWzy2ZG6TRuQMxpw5cyAqfZhc/FwD3BJXgcPMXovbLwRWx+0vSnoeOFPSS8BRZrYBQNL3gYtow3oRaYKwRo/sZtknTvFFTQviol9HcvlC0AHXs8vOb8JonCo5AThb0jKiGrlfMLP/ACYCGxP9dsZte+Pr0vbDkHQVcc3dyZOz5ZYbGoTltZtbHzfv1IlcvhDknjlr2tgmjMapgRHA0UA/sAi4P7bRl7PT2xDthzea3WVmfWbWN27cuHqNtymEBGEdNaq7Yh9n+HHRrxM33B/mj++Rti3PTuAhi3gcOAAcE7cnfQ6PB16J248v0942LM0NBAVhbb25bY8y2goX/TqQyxcI8c6845KZjR+MUys54FwASScAI4FfAY8Cl0oaJWkq0YHt42a2C9gtqT/eEVwGPDI8Q68/oWadiWN6mzAapx64Tb8OfPFHT1bs4wnUWo/58+cDfBSQpJ3AjcA9wD2xG+d7wOVmZsA2SfcDTxMd/l5rZsXiCNcA3wN6iQ5w2+YQNzS3zqJ5JzZ4JE69cNGvkVy+wHsVlvnTx492s04LsmrVKlavXr3VzPpKbi0s19/MlgHLyrRvAk5uwBCHlWgHW3kLe4enAM8UFc07ku6R9FoygEXSTZIKkrbEP+cn7nVUAEvI4a2nVnCyRqhjQrfkgp8xQmz636N8sMntZjYz/vkpdF4AS0g4+pjeniaMxHHqy/I124P6eT6d7FFR9M3s34DQMLyDASxm9iJQDGCZQBzAEttHiwEsmSWXLwRFJ950wYwmjMZx6kvhrT0V+3h6hWxSi03/ryRdBmwCbjCzN6lDAAu0dhBLmtwjC/sn+9bXyQxpq2C54GeTal02vwNMA2YCu4Db4vaaA1igNYNYluYGmPLlnwQL/lGjuv2hcDLD3BXrUgm+my2zS1UrfSXrknQAAArxSURBVDN7tXgtaSXwj/GfbRnAUk1mQQ9UcbLC0txA6vz4brbMLlWt9GMbfZFPAEXPnrYLYKlG8D0Iy8kK1RQ2dxfNbFNxpS9pFXAOcEwigOUcSTOJTDQvAX8JYGZtFcBSjeB7EJaTFaoRfHfRzD4VRd/M5pdpvnuI/m0RwBLqnZNkhDy3jpMNqhF8cBfNdsBz7wzCkofDc4dDJPjPf8PTyjqtT7WC7y6a7YGLfhly+QJvv7e/cseY6eNHu+A7mSE0n06SbskFv01w0S/DzT/eFtx3+vjRnmbByQyh+XRKcbNO++AJ18rw5jt7g/ode+RIF3wnM4Tm00nSJfiLs9ys00646JeQyxeC+nWB17l1MkVoPp1Z08a6Q0Ib4+adEr704NagfivcF9/JGCH5dMA90NodF/0Ec1es4919Byr26+3pcl/lNuDKK68EOC2ZNryIpC9IMknHJNoymzZ8aS7MG81rOLc/LvoxC1ZuCA5F/8YnT23waJxmcMUVVwA8V9ouaRIwF3g50ZbZtOFpXDR9ld/+uOgTPRSeObPzmDNnDkSR46XcDnyRQ5MCZjJteC5fCBZ8r3PbGXS86Kd5KHp7utyLoc2RdAFQMLPSwscTgaSDezE9+EQC04ZLukrSJkmbXn/99TqOujxpvHW65HVuO4WOF/2vPBR2cCvcrNPuSDoCWAJ8tdztMm2p0oY3O2X44sC5DbDiU55ErVPoaJfNpbkB3tlb+eB2RJe49eLT/KFof6YBU4En47PY44EnJJ1JxtKGL80NsCdgboObLDuNjhX9ULOOgOf/5vyK/ZzsY2YDwPji35JeAvrM7FeSHgV+KGkFcBzvpw3fL2m3pH7gMaK04d9u/ujfJ83BrefT6Tw61rwTmmrhdvfHb1vmz58P8FHgREk7JX1msL5mtg0opg3/GYenDf97osPdFxjmtOGhuXWOPqLHBb8D6diVfkiqBc+N396sWrWK1atXbzWzvnL3zWxKyd8tnzY8NLdOl+DGj3v1q06kI1f6oYEq7rPsZIk03jp+cNu5dJzoh9ryveShkzVCTZYCF/wOpiPMO7l8gZse3cZbe8KyZ3qaBSeLhGaHXdA/ucEjcVqZthf9XL7AogeeZO+B8Bzi7o/vZI3Q7LDureO0vegvX7M9leD74a2TFRas3JC6jrMLvtP2Nv3QdLIQ/cfww1snC1Qj+AvdrOMQIPqS7pH0WjL9rKTlkp6VtFXSw5LGxO1TJO2RtCX++W7iNU1PPxvqpVPEc+Q7WaEawfdVvgNhK/3vcXiq2LXAyWZ2KvCfwOLEvRfMbGb8c3WivanpZ9MkUiviZh2nXXHBd4pUFH0z+zfgjZK2fzazYkrajRyae+QwhiP97JKH063yfevrZIW0O9jp40c3aCROFqmHTf9KDg07nyopL+lfJZ0dtwWnn4XaU9AuzQ3w9nv7K3eMmTVtrK+EnEyQJq8ORIK/9vpzGjcgJ3PUJPqSlhAVobg3btoFTDaz04HriRJUHUWK9LNQWwraXL7AvSkeCi8C7WSFtILfLbngO4dRtcumpMuBPwP+KDbZYGbvAu/G15slvQCcQBPTzy5fs33wf03K4ILvZIG0gg8w/6xJlTs5HUdVK31J5wFfAi4ws3cS7eOKdUMlfYTowPaXZrYL2C2pP/bauQx4pObRl+GVFC6absd3skJo5swi7q3jDEbFlb6kVcA5wDGSdgI3EnnrjALWxp6XG2NPnTnA1yTtA/YDV5tZ8RD4GiJPoF6iM4C6p5/N5Qt0SUFZBqePH+0PhZMZQuZ0ERd8Zygqir6ZzS/TfPcgfR8EHhzkXkPTz+byBRY/NBD0cPR04bZOpy1xwXcq0TYRucvXbGfP3jCPneUXexCWA1deeSXAaSGBh/G9xXFw4XZJ8xLtDQ08XLByQ1A/F3wnhLYR/dB0C14P1ClyxRVXADxX0lw28FDSScClwAyiwMI7i+dXNDDwcGluICj61t2OnVDaQvQ9w6BTDXPmzIHI5fggQwQeXgisNrN3zexFotKIZzY68DC0jrN7oTmhZF70Q6sF+UrIqYJk4OFEIOlCUwwwTBV4mIZQs47nx3fSkOnUyqG+y54900lLmcDDwQIMgwMPJV1FZAZi8uShhTpNFk1fzDhpyOxKP02wimfPdNKQCDxcUAw8JFrBJ6OdigGGwYGHoZHmoXZ8gIljeoP6OU6RTIp+GsH3g1snDYMFHgKPApdKGiVpKtGB7eP1DjxMM7e7u8SieSdW+1FOh5I5806alMlux3eGYv78+QAfBVQp8NDMtkm6H3iayOxzrZkVfYTrEniYNtXCbRef5gsaJzWyFJF+w0FfX59t2rTp4N8fWfwTQqsfvnTLxxo0KqddkLTZzPqa/bml8zqt4LsnmjMUQ83rzJl3QgX/DrfjOxkiTW4dF3ynFjIn+iHccclM3/Y6mSI0t46bLJ1aaTvRd8F3skh3QOYGr/3g1IPMif6xR44c9J4LvpNVKuW+X9g/2QXfqQuZE/3Hlsw9TPiPPXIkL93yMRd8J7N8/aJTBq3v4DZ8p55kzmUTIuF3nHbj6xed4uLuNJzMrfQdx3Gc6nHRdxzH6SBc9B3HcToIF33HcZwOwkXfcRyng2j53DuSXgf+q45veQzwqzq+XyPwMdaHkDF+2MwGz3PcIBowr6F9/p8MJ60+PqhxXre86NcbSZuGI8FWGnyM9SELY6wnWfi+rT7GVh8f1D5GN+84juN0EC76juM4HUQniv5dwz2AAHyM9SELY6wnWfi+rT7GVh8f1DjGjrPpO47jdDKduNJ3HMfpWFz0HcdxOoi2E31J90h6TdJTibaxktZKei7+fXTi3mJJz0vaLmleE8Y3SdLPJT0jaZukz7XgGD8g6XFJT8ZjvLnVxhh/ZrekvKR/bMXx1ZNWn9fxZ7b03M7KvI4/t3Fz28za6geYA/we8FSi7ZvAl+PrLwN/G1+fBDwJjAKmAi8A3Q0e3wTg9+LrI4H/jMfRSmMU8MH4ugd4DOhvpTHGn3s98EPgH1vt/3OnzesszO2szOtGz+1hn8wN+g82peTh2A5MSEzM7fH1YmBxot8a4PebPNZHgLmtOkbgCOAJ4KxWGiNwPPAvwLmJB6Nlxteg75yZeR1/bsvO7Vad1/HnNHRut515ZxCONbNdAPHv8XH7RGBHot/OuK0pSJoCnE604mipMcbbyy3Aa8BaM2u1Md4BfBE4kGhrpfE1g5b9vq06tzMwr6HBc7tTRH8wylWjbooPq6QPAg8C15nZb4bqWqat4WM0s/1mNpNo1XGmpJOH6N7UMUr6M+A1M9sc+pIybe3sqzys37eV53Yrz2toztzuFNF/VdIEgPj3a3H7TiBZkfp44JVGD0ZSD9FDca+ZPdSKYyxiZm8B64DzWmiMs4ALJL0ErAbOlfSDFhpfs2i575uVud2i8xqaMLc7RfQfBS6Pry8nsjUW2y+VNErSVGA68HgjByJJwN3AM2a2okXHOE7SmPi6F/hj4NlWGaOZLTaz481sCnAp8H/NbGGrjK+JtNT3bfW53erzGpo0t5txaNLMH2AVsAvYS/Sv4GeADxEdjDwX/x6b6L+E6MR7O/CnTRjfbKLt11ZgS/xzfouN8VQgH4/xKeCrcXvLjDHxuefw/mFXy42vU+Z1/JktPbezNK/jz27I3PY0DI7jOB1Ep5h3HMdxHFz0HcdxOgoXfcdxnA7CRd9xHKeDcNF3HMfpIFz0HcdxOggXfcdxnA7i/wP2B/sAbeNWIwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fit,ax = plt.subplots(1,2)\n",
    "ax[0].scatter(glas_gdf_aea_rgi.index, glas_gdf_aea_rgi['h_li'])\n",
    "ax[1].scatter(glas_gdf_aea_rgi.index, glas_gdf_aea_rgi['srtm_height']);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD4CAYAAADvsV2wAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO2df5AdxXXvv+feXckoCLEWP/RjpRUyWHG05Id2DetyEkNwfpBgCMIYjPOe856NnCqSF1eSisEkKp6qcDkv5RSpelQSmeSV84Jk0A8sWxXy+GFIjMsrtKsYS0IWyAuSFv0WKyEiot29t98fd/pu397umZ6Znjtz955PlUq7d+fOdE93nz59zunTJIQAwzAMM/Mp5V0AhmEYpjmwwGcYhmkTWOAzDMO0CSzwGYZh2gQW+AzDMG1CR94FULnsssvEsmXL8i4GwzBMSzE8PHxKCHF51HWFEvjLli3D0NBQ3sVgGIZpKYjooMt1bNJhGIZpE1jgMwzDtAks8BmGYdoEFvgMwzBtAgt8hmGYNoEFPsMwTJvAAp9hmMIwfHAMj75wAMMHx/IuyoykUHH4DMO0L8MHx/CZxwYxPlnFrI4SHv/8APp6uvIu1oyCNXyGYQrB4MhpjE9WURXAxGQVgyOn8y7SjIMFPsMwhWBg+XzM6iihTEBnRwkDy+fnXaQZB5t0GIYpBH09XXj88wMYHDmNgeXz2ZyTASzwGYYpDH09XSzoM8SLSYeILiWizUT0YyLaR0QfIaL3E9GzRPR68D+3IsMwTI74suH/NYB/EUL8NICfA7APwP0AnhdCXAPg+eB3hmFSwqGLTFJSm3SI6BIAvwzgdwFACDEOYJyIbgNwQ3DZNwC8COBLaZ/HMO0Mhy4yafCh4S8HcBLA/yGifyeix4jopwBcKYQ4CgDB/1d4eBbDtDUcusikwYfA7wCwCsDfCCF+AcB/IIb5hojWENEQEQ2dPHnSQ3EYZubCoYtMGkgIke4GRAsADAohlgW//xJqAv9qADcIIY4S0UIALwohVoTdq7+/X/CJVwwTzvDBMQ5dZBogomEhRH/Udalt+EKIY0R0mIhWCCH2A7gJwKvBv88C+Grw/7a0z2IYhkMXmeT4isP/AwCPE9EsACMA/htq5qIniehzAA4BuNPTsxiGYZgEeBH4QogfAjAtJ27ycX+GYRgmPZxLh2EYpk1ggc8wMeBNT0wrw7l0GMYR3vTEtDqs4TOMI7zpiWl1WOAzjCO86YlpddikwzCO9PV0Ye0tK/H0nqO4uXchm3OYloMFPsM4MnxwDOu278X4ZBU733wbKxbMZaHPtBRs0mEYR9iGzySlKNFdrOEzjCPShj8xWWUbPuNMkaK7WOAzjCN85iqTBNPKkAU+wxQMU1ZKTlzGxKVIK0MW+AxjoEjLcKa1KdLKkAU+wxgo0jKcaX2KsjLkKB2GMcCbrJiZCGv4DGOgSMtwhvEFC3yGsVCEZTgfZ8j4hAU+w4SQp8BlxzHjGxb4DGMhb4HLjmPGN+y0ZRgLeadSYMcx4xvW8BnGQt4bZthxzPiGhBB5l6FOf3+/GBoayrsYzAwkqS3e9D12pDJFg4iGhRD9Udexhs/MeEy2eABOQluP1Mnbrs8waWCBz8x4dFv8ll2j2LprNJHQZkcq08qw05aZ8ejOTwKcnbF6HnN2pDKtDGv4zIxHd34CwJZdo5HOWJv5pgiOVPYjMElggc+0Bbot3kVo28w3ee/A3bDjENZu24OqEOxHYGLBAp9pS6SAlOYck8C0hWU2W7tWnwcAa7ftwWS1Fl03zn4EJgYs8Jm2xCXaxmS+aXaUjv681au6UalOhVKXiNiPwDjDAp9pS1yjbXTzTbOjdPTnEYDZnSWMT1RRKhHW3dbL2j3jDAt8pi2Js4tWNak0e/et/rzVq7qxelU3O2yZRPBOW6ZtcbHFp9m01exy8iTQvvBOW4aJwCXaxmTCue/Gq5sqVKPKybt/GVd44xXDBOibrIDW2GiVd1ZPpnVgDZ/xTiuaF4YPjuHT63+AiYpAZ5mwcc1H6pp1ETZahZF3Vk+mdWCBz3ilVc0LW3aNYrwSxLZXBLbsGq2XO++NVlG0wqTkm1ZUKoqAN4FPRGUAQwDeEkLcQkTvB/AEgGUA3gTwKSHEmP0OrQt3vilaLbmYbLtT5y40fE45lUfHtW8VfVLySasqFUXAp4b/hwD2Abgk+P1+AM8LIb5KRPcHv3/J4/MKAXe+RlrJvKC2XUeJ0FkmTFZEPfwxb7hvmWk1paJIeBH4RNQN4LcAPAzgj4KPbwNwQ/DzNwC8iBko8LnzNdJK5gW17SpVgbuuW4rFl140rdx5reDS9K1WWHUmLWMrKRUqRWgTXxr+IwD+FMBc5bMrhRBHAUAIcZSIrvD0rELRqp0vS4pmXrANNL3t7ljVPa3ceWrZavnKJcKRM+9h+OBY5PNbYWWQpoxZKxVZCOaitElqgU9EtwA4IYQYJqIbEnx/DYA1ALB06dK0xXHGV6O2kkbbjtgGmmz/tbesxNj5cWvb5bmCk31ry65RbB4excaXD2HLrtFIYdEKq860ZcxKqchKMBelTXxo+B8FcCsR/SaA9wG4hIj+CcBxIloYaPcLAZwwfVkIsR7AeqC209ZDeSLx3ahF02iZKWwx6q7t38wVnEkJ6evpwuDIaUxWptfBpmS0wqqzqGXMSjAXpb6pBb4Q4gEADwBAoOH/iRDid4joLwF8FsBXg/+3pX1WGtTBVJTZlvGDLiijct+Etb9+r2at4FQlpES1pGj3XF9b8ep16JozK3TCilPmvOzKRV0ZxxXMcaKoilDfLOPwvwrgSSL6HIBDAO7M8FmhyMF0YaJmC/38L16VeLZt5gDhHCrR6Ku1tbesxLrtexuEoWmg2fLcmwRpM1Zw6iRUFQJrt+3BigVzjZOOi8LiUua87cpFXBnHnSzjvL8i1NerwBdCvIhaNA6EEKcB3OTz/kkZHDmNCxNVCACTVYHHXnoD627rDbXdmnBtYB9C2OVZeQ/YvAhbrT2952hk7hvboI6z8vM90Q4sn48SEapBMsNqVTQ8XxcWPswDvNI14yqYW/H9tcVO24Hl81EuUf2UoKoQGDs/jvtuvDrWfVwa2JcQdnlWUTpcs1c9ukavCr+bexdi55tvRwpD06B2Xc6byhBXeTCVZ91tvbWjC6sCszrDy+7DPFAUu3JS8ljd5pkq2wdtIfAbBlNwDmiSxnFpYF9C2OVZRehwzV5l6O937Pz4NOG3YsHcaYLARTi4ClK1DOOTVW/ny95z/VJj2W1lTfueXeqbpVBNc+88VremZxbBLh+HGSnwTR0pzmCy4TJAfAlhl2cVwRHU7FWG6f3qwk//PY5wcBGkahmICJWqgICfzVHNtvOGPS9LoZr23nmsbouQKjstLSfwo7SCsI7kSyvy5Y2PqotLefN2BDV7lZFkktM18keeew1f/PgH63+LO1mqZeiaMwvrtu+NXX9bds4ikaVQTXvvPFa3RVhRp6WlBP6GHYewdtseVKoCszvNWoFPx1vSJWeWERKmMiUpZytvPIs7ycmBKvvFS6+fwo6R0wARJivJNEy1DElWjrbsnHHaJWsbdpYCLu298+p3ea+o09IyAn/44BjWbttTd7yOT5iFeVRHkoOka84sPPSdKc1s472NAz5rG2ESDcd23F7ccrbbxjM5UB957jW89PqpmvmlIgDUTDHjE9Faf1isf5L669k4CfHapRk27CwFnI9759Hvit7Xo2gZgT84croesgYApRIZtYKwjqQOEml7BWrL/K1K/nP5vCxthEk0HNuu0bjlLEp0jwu+VjR9PV344sc/WI/gKZeopuFPVlEF8P0DjVq/uvnJJdZfL0dUGVev6sam4dGGw8njtEuz2jBLASfvK/txUfvgTKJlBL66LJeD0bYEtnVSdZBAO7xdwC3kyjaQ4wqhJBqOrUxxJ45WsUVGrWj0HalR6O8cAB557jV8/8CpmuCsiEDnb9z85BLrH7U6BDBtRbDx3gFs3TUK2RPjtEurtGEY7bqPJE9aRuCbBGTcDqMOknK5hGq1ikoV6CwTehfNiwy5CkvElaTjxtWebJNE3ImjGbZIH/blqBWNKpTl9XF3JetafyVIlQxMbX4aWD4fHSXCREWgXKLIWP8tu0brG/0mJqvYsmsUW3eNGvvHluDzrUFStMc/3zgJ2Mrt0oZF34Wtt++WXaOFLu9MoGUEPjAlIOVh00fOvOecE0V+X9fwbDs2TSFXtmV0lstrk61Yv3eSZberYznJAPSludm0WH1Hqi5QTRuhwtImqH1i/7Fz5s1PRAAEQIQVC+aGmg03D08J7HKJQDCb3Uz9ZmD5/IZJIMx8FNaGraA9NyhgJcLm4dHETvR2ZvjgGMoXz1/gcm2hBb5J4KiROp1lQke5hErFLScKYB4k0ombdKOTy/I6aSRNXoM2zbN9TYA2LVbfkaoK1PEJ80aosDKpfaKvZ/rGrUdfOIDJSk1jr1TC468HR2qZLYGaI/bO/iVYvaobW3aNGvuNunIwKR6q+Uh1Luflp/G5aujr6cLaW1bi6T1HcVFnGc/tOx67vEVfxWSNHKfln+pa7HJ9YQW+zQ6qRupMVgTuvn7JtFOKbKYAk/MvzhZ5mwCymZvUlUSc6IuwVUeSQZBkUNie7XIvn/Zl0wStb6IDUBeoFGj/ernj+GT0Z6axra8ODlWRgu3m3oWN9VFWDqbvS/PR+MSUc3nnm2/HMl/6svH7VkCGD47VVy8dJbPy1szytCJynLoewlxYgW8T2nqkjumUIr2z29LJ6s9wya9jW0arn+sdcfWqbifBHZUnJspxbLtHVESJCZPACBtgeplcfARptDO9HaI2Qskybdk1Wh8brgIjjs/DNvnLNvjBT2r9+J7rl9ZXA9LWv3XXKB6+/dqGyUFObg3OZYujWDf9+fbTDI5MJSG0hUXHvZ8cF7UjJqcrb67fL3q0WVbIcQqhRaFYKKzAD4tI0SN1dPTObusYWUU66M8juEXSmCaguI7quBElNkwC49EXDli1fpt93EbcmHNXwRVmXwdQt/Vv2TXqPBHL9+EqTPRr1TZRHc0Dy+ejo1zrzwLApqHDWLloXn1y2Pnm2/UUyapzWe9DzUrr3DVnVt03UQ1+T4M+/kzKW1R5SsEKqVUjldIix+n1j5w54nJ9YQW+TUMxfeayLDcJ3KRaUJQAMi3rZZx1XFNImPAwCSmbSSDJpGYybXSUp6KcVMd3VnsBkqaJNq3Ukk7EaRlYbk59fN+NV+OTfd3YuONQzT9QFdYJOqyvNsteP3Z+HCUCqgIoETB2fjzV/dOsQuSqqVKt+T/W3rKy8Np9Vv6Gvp4uVN49fczl2kIJ/AMn3sWafxzCFz72AauGon/mkpMkrGPF1YJcBFDU4JTXyPvFMYVErUpM90ibNK4BuXJUVpBJVkqu39mqhTiaIrHCorVMzxyfqNn6Vy6a5zQRJ0FvV1vq4ztWdWOr4tANm6BtfbVZ9nqfz1HfT9w05cDUJFfbUiNSTz5ZUxR/Q6EE/nsTFTzz6nF8d/8JPBEI7qhZ0ZaTRCdMsOsO1rDnuWpTponJNS1CXBtx1LN9Le0HR05jMsgMWVEO6Eiiqbl8Z/jgGDYNHZ4KcSybI7FcHX59PTXnqYziWbd9r3U14Ct1r7o5zDTxyjLpNnvf7zJO2fVJVE045+M5Ps4TyMocmxVF8TcUSuBLJis1YSJjosOSpZlyksRBFxxRCbWSdjR9M47U9JPa1/PoLGF1T1KmqO/ICQaotesn+6ZsvEkdfmPnx+tRPOOWd+4zda++OezImfewZddovf6qQ1e12ft+ly6YJlF5ePpLr09FByXRyFXU92MLo43C5yTXDIoyQRVS4HeUCV1zZjklSzPlJIlC1d4aTAZKQi1bOGeSjmbajKM6oaM6QVFijZs9yExOvbC/6SsqUzm75sxC0KVQFbXf9WtNk3Ocupps9lt2jWLzcM1hDACbhw5j45qPNAi/CxPTczo1E9Mkevjt81MJ5zxppmrb2cJoXVAnf/X3IlKUCapQAv+izjJ+7WeuxBc+9gEMjtiTpekDdOO97i+yQYsJ0itMmQwIJaK6acAWzhlXm5Lhd8DUZhz5fRezhi2dQx6dp5mriyjfi0uSPF1rHDs/jiDyHSUAe46cxbrtexsOuLdNznHKbdocNhEIe6CmXMiy65E6q2NGq/hiYHnjRjA5wSZ1+ttQ2y7peQJAceziruS1MlcplMC/+oqLsf6/9td/tyVLi2v3VtFtaRIC8Kn+JbhDceLZ7G5xha0pakfiYtYw7UdIkhI5aTREnlpJ2Pux/S3MXjqwfD5md061BQENB9yv/95I/T765BwH0+awTYqG31mm+jvVI3VyjSfXNoKZJlYffUJtu6RBBUWxi7cShRL458crePCp3RCoRS+YNLi0jawK33K5BAhRS9OgmAXU+8XZfGQjzXLOZPuL+w7SHLZSpNWFK1G+Bj2f0hM7D9dNh0LI/DfC2UQIuJ1eJrNjnjh3AVfMnV3/mx6pM7B8PjbsONTgyG0GgyOnMRGsNOReBX1Fm0artr2jpJpvUezirUShBP7IyXfx+I5DAKZsnLqDKG0jmwa8baDKa9XdmUknnKSd2jZZxHkHSVcqvlYXaUmaelptN/3v6n30A+7jRo2YonJWLJhrfU/y8y1Bhsy+nsZInf3HzuHLT+0GAHzv9VMAkFjox0mzoW6sAoAnhw5P843Y+kRU+0RNFEnPOCiCXbyVKJTAVzubtHHahHCaRtYHfNg99h87hyd3HkYlcLzZ0h34wmUTWdx3YJokTQNw/7FzDVqlj9WFj/eRdIJRd9WGfc81FNImPNUQRhmVc9eHlxjfky1aqyFS58q5Dc99es/RRAJff3e/+5FleOylN6xRb3ose8UwBvU+YfNz6YT1G1lO6UcxnXHge3Xgk6KveFUKJfBVbUzaOE3EbeQ09ms9UsiU7iDOc+PkwQkTUnHegWmC0NMk/N2//gTPvHocQKNWmWZ1kVVOfJd7xY20iXqfYTmKOkoEIqpvRqsGexWkD4qI6tFAJoewXscrL3kfgLP1Z9/cu9BaJtcV2vhEFeu/N1KPUDJFvdU3pml+Bv09uaQt0QlbmQ+OnG7wo8gwVh9mpKyJW7YsJoeWTY+8/PKLcfv1S+s2fKlFyQMhwnJt2F5k2jS/pkghF2EbZ6OV+rwwLcjnqkYfgMff+c+G621aZdjqImmG0DDqQmhiSnBGoQtWIDrvS1zzlpoCoVIV+JUPXYkXfnyiHpVzx6pu9C6aV99Hsm77XtyxqtsaraW2xRc+9gHcsOKKUBu+S7/Wwx/lwS6A+YhQk59h/7FzkStOFwXAZEp99IUDGFheO12uXKK6YiVX0778dlkSp2xJZFFUv5Tp4lsyPfKcWWU8fPu19d+HD47h018fnBa77KoZDx8cwyPPvZa4s6gaT1iyNhM2W2ecPDhqGKpv56k+APcfO4dXRnfX/35z78JYSbn0a++IkSE0ytmp745VtT/bu5eCFajF24d9L67wNKVA+L2PfQC/F4QTy7pIhUGuMqTWLwMGRPBs2RbqSVf3XL+0QdDr7ynKli5/lv4IGf7okngQQN3EIlDLm5MkjYjp3rLf6u973W29+PNv7UZF1Ey7m4dH6wpekZ2zccrmO9iiwQLR6umRgamoAYnNrh/mXFQ7bdzNTa7abJwlbNw8OFH1S6MtqIJb/q9qlbbsmCb0MqrCLexQGJct9uruWJeVj3z3su2jzDquTm2XFAj6BK7HtcuInE1Dh/HNlw/VjzYEph93GGbSCLOl23aMu4Y/yvdRz4wZM41IFPr73rJrFIsvvQg3fehKPPvq8YaDZuS9w5zweWIbs2H9Mm2whfr3qltW5DqFE/iqCad30Tx0RtgUgfDQRbnB5qNXX2Y9KShK8CSxtds6Qtw8OFH1kzs0bTmE4pRZ1yrjdFD9Winc4tiYpflDd9zFXfmoGvOmocP1sFtb+V2c2qrNPlYKBENc++BILWWE6+rPNPDvu/Fquy3dsGPctjILex/y0BWTspRmdam+b/Vow44SobPDnhPJ1QmflOGDY/VJJc7mN/29RvXLNMEWpr+PT1ZbMx/++fFKgwlnVpnw0K292Hvk7DQbvkuWSfVlhR0LZxI8Ybk9XJdmpgFmm0CiJgFT/TpKhPGKmLYEjqqjq2kr7nLddG3Yd0w2Zum4+7Nv1cxL91y/dJp27WLble/ZJRumqexffmp3g9PX5VwBk9lFPxqxr6drKo+7EE6rP9vA1/uSKkSh7BiPawJR30fXnFnO5wMnuf+RM+9h48uHInMiZW3H37DjUN2kBNQ2yW28N9mk4tIvXYgaf+rf/+Cvjux3uWehBP67FyZR1kw4Y+fHG+z6QPgMKokjrHTBE5XbI61N0dW5aTLBqPW7s38JNgQ7NOXS2FbPpGW2dVCX8NGo7+iCRY2Iqgrgwad248X9J3DDiiuM2rVLnVwHmHqdKZom6lwB13TCMllaVQiUtDzuSdJI2K4BomPjTbikLU4rfPX+r573a1Nawpz3aQMapC28oujIaSYVnz6HqP4r//774+/9h8v9CiXwL57dgUkHE04aDduELnge+s7eumNNNR9ErShccHVuumhRq1d1Y9PQYSctP22Z1e8ljTZwcQD/2bd210MHBYBnXj2O54PoF5OJIouNN6rTV0bTRMXqh5ldVGeseh2hMY97WH916csmpSAOru2aRqCZnuHShnKVpzvvgfTRYIMjpxsimIBwf18UWfVLHxRK4M+ZVcbf3TsQGYaZtMOFaQJysAwfHJt2yEfa/D3q8x957rWGc0H3vHUWHSWaZmd2mdR0LV81G5hwKbOLcA87lMSGy3ekzV4V+kAtrr1cIggx3RafpB2i0PuXTK8Q9qywPrlp6DAmKgKbhw7joVt7U+1j8BHHbbuH3ue27hpNvOKwPdN0WM19N17tdA+T8x5AgzlW5u2P824GltfyK9VWD8BNH7qyfghTUuKsjJtJoQQ+4Gbjdu1wSeLCpUNNoGZSkhEEae2HquAUqGmOVQC73zqLjnIJd123pGGCs5kD9DqvXtVdi+6IEafuUkb5jkxCwHYoSdh9Xb8jhb5qT+3sKOGhT6zEniNnmxKlEUegqW1i+o5+QM+eI2ed+27cfRxhZXOxvzc4U8slbBqq5RgyPSvORKs+0/WwGhO2SVV1MH//wFTe/jhjdPWq7tjO2rgUYQNZ4QS+RHZWU/ywdORFdfYkceFdc2bVhYo0kzz0ifTpFFTBWSLgirmzceydCzVnVaWKxZdeFKpFAfbTseLGqYeVUV19yGerdZdOVWD6oSRh943zHWk+USMmgKn0A1lFaai4roZMdnupefb1dBkP6JH3Va/TMa3wAPcDc2zCJcqhKPvcW2fewzcDZ6przLiLuSvOYTU6ton48c8P4JHnXsP3D5yKrZTp70lNlhe2mTOJlp6149mFQgp8uXtMhujJ6A319KCo5ZL+cl3jwtdt39vgvKlU4qVTMJWlr6drmvZ06t0L9WttGq8qdMJi4sPi1OOgJs+qBr+bJp6tmpMtClPIZhS6wI2zJ0Anq2W03se27Bqthw5KIWs6oCfJRq+oSJ6ossn3FWUOVU2bsp3LJcKRM+/VN4npRNXH1P5J28E0Eff1dOGLH/9gorz9thQcYZsdk2rpPp25SUkt8IloCYB/BLAANTmxXgjx10T0fgBPAFgG4E0AnxJCjEXdT89fU6kKKGlKUFXyhYeZbPQkZ3HiwrX61QWfi+korFPo2hPQqPGqqxo9FC6ss4RFMESVVWXs/DhKhPoqxHQwdBL7rQ8n1sDy2kEhujM9iiyX0XqbqLn15QrpvhuvnnZAj8vkJVdueihqksgzPZTT1UkqNzttHh7FxpcPWVdWUZqrj/ZXibtJMuw+tsNubHVKo6VHlTGuYqLLPxd8aPiTAP5YCLGLiOYCGCaiZwH8LoDnhRBfJaL7AdwP4EtRNxscadw9Jk8heuylN+p5SqQ9O8xkY9PKXQcJESBAoWYSuVFDbhyx2b3VqBJde5KTkVofuVN6dmfjJjDbbkNXs46rJibNZzLhl82UNHxwrJ4PJaqDxrH5WlGc6aYcLzq6k9z3MlofwPuPnZu2QpLXhWm8YatN00avtJNsnHvIiKUwAVffVwD7Bjcv7Q+/E7gpGitKwUqrpdveQ9x66dfTrIt+yuX5qQW+EOIogKPBz+eIaB+AxQBuA3BDcNk3ALwIB4GvCx1ps//VlQtCtSSTySZuJ1MHibohxNTR5QuXwgSYuq7BfGNYDpsGo1ofAHUtUc/Tbttt6GLWUSei8cnpEQ3qxKEm/EoaNmoijf1TOtMnKyJyc5zuJC8hXaidDbWPDY6cjlwhye9EaaM+7L0+hGyUgKubQQPzq7qvIAts78W1P6r9zxaNBYRvJMwqFDhOe6vXX5ioonzRJU4d26sNn4iWAfgFADsAXBlMBhBCHCWiKyzfWQNgDQAsXbo09EVH2QWjTDYuqNqruiFE7+iq7Q+oaQjqRBO1HLbVp+7YRS2joSrE1d2eejoFF81Df8ZLr5/CjpHTuLN/ST06QU4cUiO2+T6SCKQ0SeAaV1/Rm+N0J3lYag1fxNH+ooSx671c3l0aU0GUgJPvWaC2u982yfnC9l5c+mPcPQC2NjJ9ntZPFHflIE2cdYXmfRdf5vIcbwKfiC4GsAXAF4UQ7xC5BdAJIdYDWA8A/f39AkgXmpnVcl3X7lXbX0cJuOvDSxtCulyXwyoyNGzlonkNGQ5lJ7i5dyF2jJw2brRSJxnbm5fXPPLca3jp9VO1VURFYMOOqQlpYLk54Zf+Hmwd1DVaQ408cdHM1PYIO/ha9YO4ptbwhU/tz+VettBNF3+S/L5r2KbLRNwMR6TtvbiUw9T/ovYAmPxzru8sziQQt+/09TSeheyKF4FPRJ2oCfvHhRBbg4+PE9HCQLtfCOBE0vurLzROaGZabPfXbX93fXjptPQPQDwtTQ8Nk8/Vd3fuPXI2dKNVVIKpvp6piAZTNsmB5fONCb9M99E7aNxoDVM0VdikqJbDtOt1+OAYPr3+B5ioCHQGeZjiHFXoA599MupeLlFCcU0gcVduroIqLBVo1bwAABpwSURBVCDBx3txKUfcySksgV7UOwPi7wCO23fUs5ABN7nvI0qHAPw9gH1CiL9S/vRtAJ8F8NXg/21Jn6G+0KjQTJ2kS62w75lsf2kiB8IGmd4JwjZauQ5W3eSkboKRk5k+objkzkkareErx49pk9NXDBNx1jRrN6UpSkh//zYhZ9v5nERjjxJUur8rLL++/r240WA+Jif5bP0sDVsCvbqpVBmTPnwwLvVVkqe95vIdHxr+RwH8FwC7ieiHwWdfRk3QP0lEnwNwCMCdcW6qO1dKgd0WaAzNDPsekCzPhmlm19Mlu2yKcu2wUU5elb6eRsfq2m17AKDhDNrxyejwTFkPk7nGtMPXV44VfVDGXcqGvdNT5y40/J5H7vQsw0CjzJoAjH6nO1Z1TztFzrbz2adZSqLa+YHo/Pqyrlm8Rxct2jRBmQ69Ud+ZHiWX1dnXpj7Q19PE5GlCiJdgH1s3JbmnqbHX3dZbe6FKaGbU91x31+o0RLMokTId5RI+2dfdYDcHzJuCAPfJRte4dSev3shj58cbUgmrKx49yiZqJeQigMPixuM4+FzLYCNsIgaAF/dPWQ07y9QQddEsstLsXJPPhSkicsObjHiSfGhB46HppvZI4xxWNWBbfn0d3+8xzmpBnaBKAK5dPA+9i+dhxYK5DX4kOc7lmFQDCZJs1nSpQ9pJsJA7bW3OlagTe/TvueyuNWGLCBmfrGLjjkPTTiTysaSTtkDdyQtMnzgGljeeAVoVUysePcomzCQTVpYo2zvgL6mcK6aJWIYD3vjTVzSkb1BjqptJVk7MOOa6MEVEN0HIfE6feWwwNJQxStCEXaMqAq42fJ/vMa6g1Ffc+46dw+63zjYENrikwvY9Fhr6f8JJsJAC39bYUS9Q/17SUE29g67bvjf0uDyTmeXzv3iV8+7XsHrbJr/6iieIRZfvSL9H15xZePCp3cZEWK6TgE1zb4ad0vZ+9ANTvvvjEw1ZR13SN2RBFiYRIJkADBtHcfLPuLSzi//G1V4P1ExTv3zN5bhs7uxUqRhcy6+itqG6H0eGQpuSKeonkGUxDrrmzKrv06kKJEqUWEiBn8Y04CtUU+2gKxbMjTwuTzezPPbSG1M7hLX83bZ66RON/N80aMPOVLVNVkBjNEdYNsSw9yFJqoUldWrKifXpPUexcuEleOylN+pavRACd354KRYlSMrlmyxWOUnGRNh3+nq6cHPvQvzgJ6cRtkMWiLfHI2r3sKksDdk0yyVUq1VMBhlOZqWcvIcPjuGVw2dAsG++CwtIkP4OGQr9xM7D+MTPLjTuLM5ydQvUZEwQP4cS7Bv7wiikwAfivzy10Wwn9bh8N2yZbDouTw03080se4++02DXM4XNAY0TgLyvy+HetnckP5dLenVzmDxDdGJy+g5hfRCq5fKZvyRseR3WFnq6AXVC1cNZZyJJBIrtO1E7ZJP4ZsJSDIe1u8kUK0mzcpRhujJyS9YTQD0lCBDua+vrqZ058fiOWu6rSlXgWz88AqB2xKh8b82IzBpYXsvbn8bMVViBH4c0zow439UHj/5dXQDpXn09bM42Aej5X8bOj8eexIDp+c0/2VcbkBtfPtQwCZg2L7nEH5veSRRhy+uottC/O/eiTjzxhY80JQRypqE6JvUdsnF9M6Z9JLbnmdpd76eqhp/Gfj84cjo40L1GtSqw98jZhr7sEtixelU3nth5uMHRDeW9+XCm2jZ0hUVlJenvM0Lgp7El+/zuOxcm8akPL2nQctS87isXzQuNm/67f/0JvvvjE3XTkLoETaJBmDqImjKirEUdqfVSsz66HODtStjSP6otkjjGmhUT32qEhQLHHRMu14e1u95PATSchZCmr5XLhMlA6AugtlM9ZmBHX09X3WdmGptpfVm2HdMuUVlxmRECP41H39d3palEZs1UtZzNSj70hz7RGEqoCt/n9x2v5+InAB+9ppb/BUh+bqfeQVy0BD0v/sqFl0yLP05jhzc9f/jgGN46817oaUhxNZwsY+JbHfkuTaHAcceE6fq42qmpn/pg5cJL8Mro2frvB07WwtVlaKhrYIfqMzNFGqWJKDJNGED0YTd5pUfOnTRLHV/ftWXXlCYboNaAe4+cbUjDIL//1pn3sDGwEwI1e6PM/5Lm8A9bucO+r+fFn3tRp9NGs6TPVw+86SwTbvrQlbh87uxEZVdpdhRRUYgTfWUKBY4bcWLS0LPQTsMw+ZxkGUxcu3ge1n5iym/hUq4wn1kaU4tpwtx/7Fxoyunc0iMXBdfO5JIiIMlzbdk19QQX+u/q97cGKRNKpVq+IJONsxkJqmymEwD1ycmHIB0+OIa//def4Ll9x+up7icqom7WemLn4XreJF/1mOmkiTlPE3GiXu9bQYnCVGd1slePFJX0Lp7ntUxpZYg+YUalnNZP6irNmjPXcOtpzBiB70KWS3zbLH/Hqm5sHjpcT+hlCzEL0xLU5bdruoA0tmub7V8NnVNj3pMI0g07DjUcVC4hgnUXcVx8OLlajTQx51HvKGmqkFcOn8GDT+3OLIrKVGd9IvsfN30QD31nKsPqykXznA/v8U2U0qlG15lSTuvZesslQnX8/DmXZ7eVwM96iW+a5ft6urBxjVskSZSWEJUJUxJnYrMNYr0s6rurVKq4+7rkMe/DB2vHWE4T9gBu/blF2P6jo8ZdxGH3c7UNz3SSrGpc3lHcaDapoDw5dBjPvHocALBpeBQb7/XvR7GtSHWtWUaprVw0zxp1ljUu7zGqDaUZDpjaVf6VZuXSaSXyWuLHFTomARZnsnK9Ns4g1t9dGm1tcOQ0KlXduFXjX/Yeq4e3VqqiftSiDXbM+slnFKW9q31qfHL6aWk60n81qczqWZp39ARxsgz66nRWR6kWfdYkk5P+Xl3GZlQbmsbiVxzL01YCvxWW+DYBFmeycr02ziTi893JDSS1tBPAVZdfjJGT7zbE16upI8KSwLWbYzbMOWmLmY8S5nG0TvmuX3r9FHa++XZo/p1NQ4cbPstCyTIlTNT/rqc5JpgzwvqeJE3v1XVshimJacZiWwl8wO8SP4sYb5sAi9PItmv18sZd8fh6d6bl9mceG5yWQyjqGEOgvRyzUc5J247pKGGua+9hWqd6WlpU/h11o9LPdTdGxbjWN6q/h9XfluZ49aruhl3zQPyos7jv1XfOnaRjse0EflJcNau0RG1QcX2GScszlTerFY9rugpJ0oNRTJNHXs64rHFxTprsvVErINekXH09U6elReXMOXLmvQbHfhJhn/YMBll3uVHKdK6xDKmOu0pMs9ksT78SC3wHkmhWSclKCIetHFwduq4rmiSToV6OuCsaWb5Pf31qpZCFgzBPXJyTcR2AQLykXGHPk6HJcgNiR7mEu65bkijbpev4CiuPXndV2DdEnZUodMOfCZf3WkQTMgt8B5JoVmnIQgOIW944+XR0fE2Gcd+Dvslt667RQgwyX9gESFp7r/SpxDHtAbV2lr/r5hOgFs21+NKLErVBnP5qq39Y3RuizqoCd123BID7KWmuwtxUtjzTfbDAdyCJZpUW350ibnl1oR0nn47pfTWjk0dtcpsJJFEGor4Tt2+ErXjDkvLFLXMWtm41u63aR3uVUM2osGfbvV1Im+gx7ftggR8QFcsdV7NKW5YkTqQk2obtu7rQtp3naXuOz1QMrrhucmOmY9LabUSteMslwp39S1JvtPI9vsKOyGxWtFfS5/jyGbLAh9vL9Nn54sQ8u3SKLNJDmya5qCMmTQmzgOZtte/rcd/kxjQic8fLyXLjmo84reBkpk3A7HgvEvq40tOONyPaK6kpeHCkMYNt0jHEAh/NjeW2pUJNEy6Zpvxh3zU5Um33lQnQ5FkAar1sp3aZ3k1e4WrtzpZdo/WDQsYrAltC/B9SGfjbIJ23nmnTZZXgSpo+ERWG3DVnVkNEVzMmrKTP0TPYJjneEGCBD6C5sdy6gDUdghK3U9jK7zJY0tRdtYeu3banHnc9bqiX7dQu9V7tvmM2T3RnpYvz8oUgyR0w5SSXSb3KQQLApInvgGxWrnJcnXtvYpqC0ixlIclz9Ay2SY43BFjgA2hu+JQuYPVDUFzCJV3K7zpYktZdvX8pOFBc5dS5C6HLZ5122zFbNFav6sYm5dwG06lVKoMjjekxSkQQQN3skDbxnXxGUjNGWBjy/mPn8LVn9tf3Htg2mxUJX0opC/yAZs7uptN90jakXv44AjRJ3dX7CyFAhHqKYyGAF/efiBXb3E47ZotIX08XNt4bb1Up02PIdN4rFszFk8pRgC6J78JIY8YIW/XWNPupa0tEhe9vvpRSFvgx8BVaqAvYLFYXWQtQef/xiSqqABAsNYWohUPK2OaT5y7gxDv/if3HznkNDWT8k3ZVCaAhB9KslP0ujRnDVr7pKxM0nD3RbOLIFB9KKQt8R7LOpe+7w8UVoHEnM3n/dd/Zi1dGz0IAIFHLzS1EbTv9JbM7sCE4xeuV0d0AEGrTZYdra2FqL3kUYJyzG2zEUVpcDzYyrUzS+BnSYJMpWe5ZYYHvSCvYmE1hkVnH+O47+k79546OxjN7H3nutYZrn95zNLfBxWSDTTi5nt0Qhqq0dM2ZZY3+SZKrvwgrSZNMAbLds8IC35E0kTDNII3QTjqZDY5MZUQk1A6YuOf6pfV3snLhJfje66fq169ceEmiujHFxNbnfCpH8nthfTvu84qykjTJlKwVSxb4jqSJhGkGaTpKUnu//r07VnVPeye//fOL8O1XjkAI4B++/wbOXZjM7Kg7prnY+pxv/9FW7fxWvW+3qsPfttrIsi4s8GPgEgkjP2+2xp+m0ydd5pq+p++qPf0fNUebQG1Dz4Ydh1It85niYOtzPs0m8iCV+vmt5el9u0hmGhNxjuDMui4s8FNg2rmXl8aftqMkXebq37Pl4JEaWtTBGUzrENbnfJlNBkdOY0I5JvGTfebVYVHMNDo+UoX7hAV+CvQOn7djtwid3iQEZNTG5uHRWDnHmeIT1ud8+LfUWHygtqFv+OBY7v3clbxlgg4L/JToHb4VbYm+MS1T+3q6cIdyrFyrDNiZQrODC3z5t9RYfAB49tXj+LfXT7aMSbBo/gUW+B4pui0RyDeqSN38ov7OZEsewQW+NFspMFvVJOgqE+RpYQRkGtTAAt8zRTCr2Mg7qijv57creZgVfGm2UmBu3TWKTUOHMVkRIKLE2SLzIEombNhxCH/+rd2QropNw6OZHc9Z8n5HDSL6DSLaT0QHiOj+rJ/HNDJ8cAyPvnCgrtmbooqahfr8CxO1TTlM9kjhW6Z0p1DFQQrqP/q1Fakn9r6eLjx8+7V46NZelIKD0ddu21Pfxd3KyNw+il+6vmEtCzLV8ImoDOBRAL8KYBTATiL6thDi1Syf286oJhsA0074ydOeOLB8PjpKhPGKgADwxM7D6F00j3ffZkxepkbfq92x8+OoVIW3bJy+SWIu1XP7SDYPjyY6/D2KrE061wE4IIQYAQAi+iaA2wCwwM8A3WRyx6ruaSmK8/Qx9PV04c7+JXg80MwqBRy0M5UimxpdGVg+H+US1Xd3V6rhB7U0k6TmSjW3j/RRALUD4FXTmy/fW9YmncUADiu/jwaf1SGiNUQ0RERDJ0+ezLg4xUE1tZh+T4JushHAtKV8X08X7rvx6twGyepV3egoTaXVkil0GSaKvp4urLutF+Wg+wjUNOE0Y8YXSc2lcvX1x7++Ag/ffi3e1znd9CYnk689sx+feWwwVX2z1vBNCfMa1i9CiPUA1gNAf3//9LXNDMR0mPK67XtTOzNNqQ6KFgopB62vFLpMe3HP9Uux98hZbNhxqJaGu1KMiJ20O91l+U3nRvt0umct8EcBLFF+7wZwJONnFh69AZ/ec9RLg9pstXkPBh2ZQrdIExHTOqxe1e3l0CCf+PKTmExvPmP5sxb4OwFcQ0RXAXgLwN0A7sn4mYXHln7AR4O2iq22VcrJFI8w4Zr3PhPXZ27YcQhP7zmKm3sXRgYt+HS6ZyrwhRCTRPT7AP4fgDKAfxBC7M3yma2ALf0Aa7xMq9MsgWsSrq2yz2PDjkP48lO1A4Fk+nAXoe+jLplvvBJC/DOAf876Oa2GLf3ATKQoZwYw2ZK3wC1a3hobT+85Ou13m8D3PXZ4py2TKXkLAaZ55C1w1XOWAeCVw2cKmWjt5t6FDQcD3dy70HhdFmMn8522THuT9+5epnnksaNXpa+nC2tvWQkioCKAZ149jk9/PV0YYxbcc/1SfOX2a/FL11yGr9x+rVW7z2LssIafAWzCmKJo2QKZ7ChC8sCx8+NQN64W1bRzz/VLI+32WYwdFvieKaIJI+/IhbyFANM88vZFDSyfj3KZMBkkp+koU8sqGVmMHRb4nsnbjqlThAkobyHAtBliSsVv9Z2cvscO2/A9k7cdU4dt6Ew7UUtGNvV7pcKpO1RYw/dM0UwYJjsg+xiYmcrA8vnoLNcysgLFULqKBAlRnEVPf3+/GBoayrsYM46wlMlF8DEwjE+adXpUkSCiYSFEf9R1rOG3Aaod8NEXDhTKx8AwvplpPiOfK3IW+G0Gh0kyTOvgO+iCBX6bUTQfA8MwdnxH/bHAb0Nm2pKXYWYqvlfkLPAZhmEKimlFnsamzwKfYRimwKgr8rQ2fd54xTAzDB/nIzPJyPrdp91IyRo+w8wgipBKo11pxrtPa9Nngc8wM4ii5XJqJ+K++yS2+LRRdizwGWYGwfss8iPOu0+zGkgTZccCn2FmELzPIj/ivPu8VmIs8BlmhsH7LPLD9d3ntRJjgc8wDNNk8lqJscBnGIbJgTxWYhyHzzAM0yawwGcYhmkTWOAzDMO0CSzwGYZh2gQW+AzDMG0CC3yGYZg2gQU+wzBMm8ACn2EYxiNFTk/NG68YhmE8UfT01KzhMwzDeCLuASXNXg2whs8wDOOJZqVITgoLfIZhGE8UPUUyC3yGYRiPxE2RPD5RBRGha86szMuWyoZPRH9JRD8moh8R0VNEdKnytweI6AAR7SeiX09fVIZhmJlDX08X1t6yEqUSoSoE1m3fiw07DmVq00+r4T8L4AEhxCQR/QWABwB8iYh+BsDdAFYCWATgOSL6oBCikvJ5DMMwM4ax8+OoCoGqAMYnqli7bQ+qQmRm00+l4QshnhFCTAa/DgLoDn6+DcA3hRAXhBBvADgA4Lo0z2IYhplpSLNOmVDX9F0jfJLg04b/3wE8Efy8GLUJQDIafDYNIloDYA0ALF261GNxGIZhio3q5O2aMwvrtu/N9NjDSIFPRM8BWGD404NCiG3BNQ8CmATwuPya4Xphur8QYj2A9QDQ399vvIZhGGamojp5VyyYm+mxh5ECXwjx8bC/E9FnAdwC4CYhhBTYowCWKJd1AziStJAMwzDtQNbHHqaN0vkNAF8CcKsQ4rzyp28DuJuIZhPRVQCuAfBymmcxDMMw6Uhrw//fAGYDeJaIAGBQCPF7Qoi9RPQkgFdRM/XcxxE6DMMw+ZJK4Ashrg7528MAHk5zf4ZhGMYfnDyNYRimTWCBzzAM0yawwGcYhmkTWOAzDMO0CSzwGYZh2gQW+AzDMG0CC3yGYZgC4/MYRD4AhWEYpqD4PgaRNXyGYZiCEvdQ9ChY4DMMwxQUNV++j5TJbNJhGIYpKHEORXeBBT7DMEyB8ZkymU06DMMwbQILfIZhmDaBBT7DMEybwAKfYRimTWCBzzAM0yawwGcYhmkTSAiRdxnqENFJAAdzLsZlAE7lXIa0zIQ6ADOjHlyH4jAT6mGrQ48Q4vKoLxdK4BcBIhoSQvTnXY40zIQ6ADOjHlyH4jAT6pG2DmzSYRiGaRNY4DMMw7QJLPCnsz7vAnhgJtQBmBn14DoUh5lQj1R1YBs+wzBMm8AaPsMwTJvAAp9hGKZNaCuBT0T/QEQniGiP8tn7iehZIno9+L9L+dsDRHSAiPYT0a/nU+rpWOrxEBG9RUQ/DP79pvK3wtWDiJYQ0QtEtI+I9hLRHwaft0x7hNSh1drifUT0MhG9EtTjfwaft1Jb2OrQUm0BAERUJqJ/J6Ltwe/+2kEI0Tb/APwygFUA9iif/S8A9wc/3w/gL4KffwbAKwBmA7gKwE8AlPOuQ0g9HgLwJ4ZrC1kPAAsBrAp+ngvgtaCsLdMeIXVotbYgABcHP3cC2AFgoMXawlaHlmqLoGx/BGADgO3B797aoa00fCHEvwF4W/v4NgDfCH7+BoDfVj7/phDighDiDQAHAFzXlIJGYKmHjULWQwhxVAixK/j5HIB9ABajhdojpA42ClcHABA13g1+7Qz+CbRWW9jqYKNwdQAAIuoG8FsAHlM+9tYObSXwLVwphDgK1AYwgCuCzxcDOKxcN4rwwVwEfp+IfhSYfOSyr/D1IKJlAH4BNa2sJdtDqwPQYm0RmBF+COAEgGeFEC3XFpY6AK3VFo8A+FMAVeUzb+3AAt8OGT4rcgzr3wD4AICfB3AUwNeCzwtdDyK6GMAWAF8UQrwTdqnhs0LUw1CHlmsLIURFCPHzALoBXEdEvSGXF7Ieljq0TFsQ0S0ATgghhl2/YvgstA4s8IHjRLQQAIL/TwSfjwJYolzXDeBIk8vmjBDieNDhqwC+jqmlXWHrQUSdqAnKx4UQW4OPW6o9THVoxbaQCCHOAHgRwG+gxdpCotahxdriowBuJaI3AXwTwK8Q0T/BYzuwwAe+DeCzwc+fBbBN+fxuIppNRFcBuAbAyzmUzwnZIQJuByAjeApZDyIiAH8PYJ8Q4q+UP7VMe9jq0IJtcTkRXRr8fBGAjwP4MVqrLYx1aKW2EEI8IIToFkIsA3A3gO8KIX4HPtshb490M/8B2Ijasm4CtdnxcwDmA3gewOvB/+9Xrn8QNc/3fgA3513+iHr8XwC7Afwo6AgLi1wPAL+I2vLzRwB+GPz7zVZqj5A6tFpb/CyAfw/KuwfA2uDzVmoLWx1aqi2Ust2AqSgdb+3AqRUYhmHaBDbpMAzDtAks8BmGYdoEFvgMwzBtAgt8hmGYNoEFPsMwTJvAAp9hGKZNYIHPMAzTJvx/1azNnQQLu0IAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "glas_gdf_aea_rgi['diff'].plot(marker='.', linestyle='none');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "21.185592375756876"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "glas_gdf_aea_rgi['diff'].mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### geoid offset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
